Create a table summarizing data layers within each MTBS fire. Layers were assembled and reduced in Google Earth Engine, then the summary file exported to Google Drive. Downloaded the summary table locally from Google Drive, then ran this script to visualize and explore the data.

# how many fires in this data set? 
message("# of fires in dataset: ")

## # of fires in dataset:

message(nrow(fire_df))

## 37

# create output directory for figures
if(!dir.exists(here::here("figures"))){
  dir.create(here::here("figures"))
}

# create a table to summarize attributes of interest per fire 
fire_table <- as.data.frame(fire_df) %>% 
  # select columns of interest
  dplyr::select(Fire_Name, Year, Acres, gedi_coverage, lodgepole, ponderosa, spruceFir, disturbed_burned, disturbed_unspecific, 
                #disturbed_logged, regenerating_disturbed, regenerating_harvested
                ) %>%
  # reorder the rows based on multiple variables 
  dplyr::arrange(desc(gedi_coverage), desc(lodgepole), desc(ponderosa), desc(spruceFir), desc(disturbed_burned)) %>%
  # rename the columns for interpretation
  dplyr::rename("Fire Name" = Fire_Name, 
                "Year" = Year,
                "Acres" = Acres,
                "GEDI %" = gedi_coverage,
                "Lodgepole %" = lodgepole,
                "Ponderosa %" = ponderosa,
                "Spruce Fir %" = spruceFir,
                "Disturbed Burned %" = disturbed_burned,
                "Disturbed Unspecified %" = disturbed_unspecific,
                #"Disturbed Logged %" = disturbed_logged,
                #"Regenerating Disturbed %" = regenerating_disturbed,
                #"Regenerating Harvested %" = regenerating_harvested
                ) 

# add a categorical column with dominant forest type 
temp <- fire_df %>% 
  dplyr::select(Fire_Name, lodgepole, ponderosa, spruceFir) %>% 
  tidyr::gather(major_forest_type, forest_percent, lodgepole:spruceFir) %>% 
  dplyr::group_by(Fire_Name) %>% 
  dplyr::slice(which.max(forest_percent)) %>% 
  dplyr::select(Fire_Name, major_forest_type) %>% 
  st_set_geometry(NULL)
fire_df <- fire_df %>% 
  dplyr::left_join(temp)
remove(temp)



kableExtra::kable(fire_table) %>%
  kableExtra::kable_styling(bootstrap_options = "striped", 
                            full_width = F, 
                            position = "left")

Fire Name	Year	Acres	GEDI %	Lodgepole %	Ponderosa %	Spruce Fir %
MAD CREEK	2001	1435	99.7	1.1	0.0	74.6
BIG ELK	2002	4347	98.5	19.3	47.9	0.0
CRYSTALFIRE	2011	2855	95.3	2.7	54.2	0.0
LOWER NORTH FORK FIRE	2012	3436	93.8	0.2	31.2	0.0
#6	1989	1547	93.2	0.0	6.7	0.0
SNAKING	2002	2080	91.4	0.2	65.4	2.7
FOURMILE CANYON	2010	5865	91.1	1.2	67.0	0.0
BIG RED	2017	2868	90.0	80.6	0.0	0.9
UTE CREEK	1994	3225	89.2	2.1	0.0	36.0
PICNIC ROCK	2004	9014	88.2	0.0	15.9	0.0
MT. ZIRKEL COMPLEX (HINMAN)	2002	17627	86.4	10.1	0.1	23.1
BUFFALO CREEK	1996	9792	85.5	0.0	9.4	0.0
OLD MAN	2001	1169	83.4	0.0	0.0	0.0
MEADOW CREEK	2010	1432	83.3	3.8	0.1	8.9
MT. ZIRKEL COMPLEX (BURN RIDGE)	2002	15478	71.8	4.4	0.0	58.0
BOBCAT	2000	9131	69.8	1.4	75.9	0.0
BIG FISH	2002	17273	68.3	3.2	0.0	58.2
HAYDEN PASS	2016	17978	61.7	2.0	21.5	15.0
BEAVER CREEK	2016	44222	60.4	47.7	0.0	3.6
WALDO CANYON	2012	20112	58.6	0.3	47.2	0.2
LOST SOLAR	2016	5142	56.3	3.1	0.0	52.9
SUNNYSIDE	1989	1416	47.8	0.0	27.5	0.0
CANYON	1988	1092	47.5	0.0	40.8	0.0
HIGH PARK	2012	90769	46.8	13.9	49.6	1.1
LOST LAKES	2002	5888	44.9	0.2	0.0	51.8
UNNAMED	1989	1305	43.6	0.0	43.7	0.0
DUCKETT	2011	4705	41.2	0.4	17.4	0.0
COAL SEAM	2002	12004	37.3	0.1	0.0	0.5
OVERLAND	2003	3232	35.6	1.8	79.6	0.0
HAYMAN	2002	129417	25.4	0.9	57.4	0.6
GREEN CREEK	2002	4342	21.6	16.4	0.0	55.5
SPRINGER	2012	1666	14.7	0.0	83.9	0.2
RED FEATHER	2017	1582	13.0	0.4	79.5	0.0
GRACE CREEK	1988	1455	11.7	32.1	0.5	3.4
SPRING CREEK	2002	11668	6.5	0.7	0.0	45.5
HIGH MEADOWS	2000	9607	0.0	0.0	75.0	0.0
SCHOONOVER	2002	2821	0.0	0.0	73.2	0.0

# count number of fires with each dominant forest type 
kableExtra::kable(table(fire_df$major_forest_type),
                  col.names = c("Dominant forest type", "# of fire events")) %>%
  kableExtra::kable_styling(bootstrap_options = "striped", 
                            full_width = F, 
                            position = "left")

Dominant forest type	# of fire events
lodgepole	4
ponderosa	22
spruceFir	11

Explore the fire data

Map

library(leaflet)
library(rgdal)
library(geojsonio)

# read geojson file for mapping with leaflet
fire_gjson <- rgdal::readOGR(fire_filename)

## OGR data source with driver: GeoJSON 
## Source: "C:\Users\tyler\Documents\Field Planning\CAREER\data\fires_stats_20210414.geojson", layer: "fires_stats_20210414"
## with 37 features
## It has 14 fields

# add a column with fire perimeter centroid point locations for map labels
fire_df$centroid_lon <- NA
fire_df$centroid_lat <- NA
fire_df$centroid_label <- NA
for (i in 1:nrow(fire_df)){
  # create popup label using fire name and year
  content <- paste(sep = "<br/>",
    fire_df$Fire_Name[i], fire_df$Year[i])
  #print(content)
  fire_df$centroid_label[i] <- content
  
  # get lon, lat coordinates of geometry centroid  
  fire_df$centroid_lon[i] <- st_centroid(fire_df$geometry[i])[[1]][1]
  fire_df$centroid_lat[i] <- st_centroid(fire_df$geometry[i])[[1]][2]
}


# color palette
pal <- colorNumeric(
  palette = "viridis",
  domain = fire_gjson$Year
)

# map the fire polygons
map <- leaflet(fire_gjson) %>%
  addTiles() %>%
  addPolygons(stroke = FALSE, smoothFactor = 0.2, fillOpacity = 0.7,
    color = ~pal(Year), group = "MTBS fires"
  ) %>%
  addLegend("bottomright", pal = pal, values = ~Year,
    title = "Fire Year",
    # get rid of the comma in years (1,990 becomes 1990) for legend entries
    labFormat = labelFormat(big.mark = ''),
    opacity = 1
  ) %>% 
  # set base map type
  addProviderTiles(providers$CartoDB.Positron, group = "Carto") %>%
  addProviderTiles(providers$Stamen.Terrain, group = "Stamen")

# Icons 
icons <- makeAwesomeIcon(
  icon = "fire",
  iconColor = "gray",
  library = "fa",
  markerColor = "gray"
)


# Add Fire Name popup labels to each polygon centroid
  map <- map %>%
    #addMarkers(lng = fire_df$centroid_lon,
    addAwesomeMarkers(lng = fire_df$centroid_lon,
               lat = fire_df$centroid_lat,
               icon=icons,
               popup = fire_df$centroid_label) %>% 
    addLayersControl(baseGroups = c("Carto", "Stamen"), 
                   overlayGroups = c("MTBS fires"))  


# display map
map

Fire timeline

# calculate minimum and maximum fire years
year_min <- min(fire_df$Year)
year_max <- max(fire_df$Year)

# set forest type colors for the figures
                       # lodgepole  ponderosa  spruceFir
#forest_type_colors <- c("#005a32", "#74c476", "#e5f5e0") 
forest_type_colors <- c("#04663B", "#41ab5d", "#c7e9c0")


# reorder the rows based on year
fire_df <- fire_df %>% 
  dplyr::arrange(desc(Year)) 

# calculate point size - scaled relatively by fire size
tmp <- log(fire_df$Acres / max(fire_df$Acres))
tmp <- tmp + abs(min(tmp)) + 1.25

# create timeline figure 
ggplot(fire_df) + 
  geom_point(aes(x = reorder(Fire_Name, Year), y = Year, 
                 # color points based on dominant forest type
                 fill = major_forest_type), 
             # add a black outline around each point
             color = "black", shape=21, 
             # size of each point, thickness of outline
             size = tmp, 
             stroke = 0.3) + 
  scale_y_continuous(breaks = seq(year_min, year_max, by = 4)) + 
  labs(title = "Fire Event Timeline", y = "Year", x = "Fire Name\n",
       fill = "Dominant forest type") + 
  # Put fire name on Y axis, year on X axis
  coord_flip() + 
  # set the point fill colors using hex codes 
  scale_fill_manual(values = forest_type_colors,
                    labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) + 
  # set the point size in legend 
  guides(fill = guide_legend(override.aes = list(size=4))) + 
  # clean theme for plot with white background
  theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_timeline.png")))),
       width = 8, height = 6)

Fire year histogram

hist_year_title <- paste("Histogram: Fire Years,", year_min, "-", year_max)

# create histogram, fire count per year.
# color the bars based on the dominant forest type of each fire event. 
ggplot(fire_df, aes(x=Year, fill = major_forest_type)) + 
  geom_histogram(binwidth=1, color="black", 
                 # set the outline thickness and bar transparency. 
                 size = 0.2, alpha=0.9) + 
  labs(title = hist_year_title, y = "Count", fill = "Dominant forest type") + 
  # x axis label years from min to max year in increments of 4
  scale_x_continuous(breaks = seq(year_min, year_max, by = 4)) + 
  # set color scale of the dominant forest types
  scale_fill_manual(values= forest_type_colors,
                    labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) + 
  theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_hist_years.png")))),
       width = 7, height = 5)

Fire size histogram

# color the bars based on the dominant forest type of each fire event. 
           # convert from Acres to Hectares? * 0.404686
ggplot(fire_df, aes(x = Acres * 0.404686)) + 
  geom_histogram(aes(fill = major_forest_type), color="black", 
                 # set the outline thickness and bar transparency. 
                 size = 0.2, alpha=0.9) + 
  labs(title = "Histogram: Fire Size", y = "Count", x = "Size [Hectares]",
       fill = "Dominant forest type") + 
  # set color scale of the dominant forest types
  scale_fill_manual(values = forest_type_colors,
                    labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) + 
  theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_hist_size.png")))),
       width = 7, height = 5)

Fire size as a function of year

ggplot(fire_df, aes(x=Year, y=Acres * 0.404686 , 
                    log="y", label = Fire_Name)) + 
  geom_point(aes(# color points based on dominant forest type
                 fill = major_forest_type), 
             # add a black outline around each point
             color = "black", shape=21, 
             # size of each point, thickness of outline
             size = 4, stroke = 0.5) + 
  labs(title = "Fire size vs year", 
       y = "Hectares [log transformed]", 
       fill = "Dominant forest type") + 
  # log transform the y axis to space out the points more 
  scale_y_continuous(trans='log2') + 
  scale_x_continuous(breaks = seq(year_min, year_max, by = 4)) + 
  scale_fill_manual(values = forest_type_colors,
                    labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) + 
  # add Fire Name labels
  geom_label_repel(aes(label = Fire_Name),
                   box.padding   = 0.35, 
                   point.padding = 0.5,
                   segment.color = 'grey50') + 
  
  theme_bw() + 
  # font size
  theme(axis.text.y = element_text(size=16),
        axis.title.y = element_text(size=16),
        axis.text.x = element_text(size=16),
        axis.title.x = element_text(size=16),
        # move the legend to the top of the figure
        legend.position="top")

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_size_vs_year_scatterplot.png")))),
       width = 7, height = 9)

VictoriaFireSummaryScript

Victoria Scholl & Tyler McIntosh

2022-04-26

Explore the fire data

Map

Fire timeline

Fire year histogram

Fire size histogram

Fire size as a function of year